! COPYRIGHT (c) 2011 Science and Technology Facilities Council
!
! Author: Jonathan Hogg
!
! History: See ChangeLog
!

module hsl_mc64_double_ciface
   use iso_c_binding
   use hsl_mc64_double, only:                &
      f_mc64_control       => mc64_control,  &
      f_mc64_info          => mc64_info,     &
      f_mc64_matching      => mc64_matching
   implicit none

   type, bind(C) :: mc64_control
      integer(C_INT) :: f_arrays
      integer(C_INT) :: lp
      integer(C_INT) :: wp
      integer(C_INT) :: sp
      integer(C_INT) :: ldiag
      integer(C_INT) :: checking
   end type mc64_control

   type, bind(C) :: mc64_info
      integer(C_INT) :: flag
      integer(C_INT) :: more
      integer(C_INT) :: strucrank
      integer(C_INT) :: stat
   end type mc64_info

contains
   subroutine copy_control_in(ccontrol, fcontrol, f_arrays)
      type(mc64_control), intent(in) :: ccontrol
      type(f_mc64_control), intent(out) :: fcontrol
      logical, intent(out) :: f_arrays

      f_arrays          = (ccontrol%f_arrays.ne.0)
      fcontrol%lp       = ccontrol%lp
      fcontrol%wp       = ccontrol%wp
      fcontrol%sp       = ccontrol%sp
      fcontrol%ldiag    = ccontrol%ldiag
      fcontrol%checking = ccontrol%checking
   end subroutine copy_control_in

   subroutine copy_info_out(finfo, cinfo)
      type(f_mc64_info), intent(in) :: finfo
      type(mc64_info), intent(out) :: cinfo

      cinfo%flag        = finfo%flag
      cinfo%more        = finfo%more
      cinfo%strucrank   = finfo%strucrank
      cinfo%stat        = finfo%stat
   end subroutine copy_info_out
end module hsl_mc64_double_ciface

subroutine mc64_default_control_d(control) bind(C)
   use hsl_mc64_double_ciface
   implicit none

   type(mc64_control), intent(out) :: control

   type(f_mc64_control) :: fcontrol

   control%f_arrays  = 0 ! false
   control%lp        = fcontrol%lp
   control%wp        = fcontrol%wp
   control%sp        = fcontrol%sp
   control%ldiag     = fcontrol%ldiag
   control%checking  = fcontrol%checking
end subroutine mc64_default_control_d

subroutine mc64_matching_d(job,matrix_type,m,n,cptr,crow,cval,ccontrol,cinfo,&
      cperm,cscale) bind(C)
   use hsl_mc64_double_ciface
   implicit none

   integer(C_INT), value, intent(in) :: job
   integer(C_INT), value, intent(in) :: matrix_type
   integer(C_INT), value, intent(in) :: m
   integer(C_INT), value, intent(in) :: n
   type(C_PTR), value :: cptr
   type(C_PTR), value :: crow
   type(C_PTR), value :: cval
   type(mc64_control), intent(in) :: ccontrol
   type(mc64_info), intent(out) :: cinfo
   type(C_PTR), value :: cperm
   type(C_PTR), value :: cscale

   integer :: i
   integer(C_INT), dimension(:), allocatable, target :: &
      fptr_alloc, frow_alloc, fperm_alloc
   integer(C_INT), dimension(:), pointer :: fptr, frow, fperm
   real(C_DOUBLE), dimension(:), pointer :: fval, fscale
   type(f_mc64_control) :: fcontrol
   type(f_mc64_info) :: finfo
   logical :: f_arrays

   ! Copy data in and associate pointers correctly
   call copy_control_in(ccontrol, fcontrol, f_arrays)
   call C_F_POINTER(cptr, fptr, shape = (/ n+1 /))
   if(.not.f_arrays) then
      allocate(fptr_alloc(n+1))
      fptr_alloc(:) = fptr(:) + 1
      fptr => fptr_alloc
   endif
   call C_F_POINTER(crow, frow, shape = (/ fptr(n+1)-1 /))
   if(.not.f_arrays) then
      allocate(frow_alloc(fptr(n+1)-1))
      frow_alloc(:) = frow(:) + 1
      frow => frow_alloc
   endif
   call C_F_POINTER(cval, fval, shape = (/ fptr(n+1)-1 /))
   if(.not.f_arrays) then
      allocate(fperm_alloc(m+n))
      fperm => fperm_alloc
   else
      call C_F_POINTER(cperm, fperm, shape = (/ m+n /))
   endif
   if(C_ASSOCIATED(cscale)) then
      call C_F_POINTER(cscale, fscale, shape = (/ m+n /))
   else
      nullify(fscale)
   endif

   ! Call the Fortran routine
   if(associated(fscale)) then
      call f_mc64_matching(job, matrix_type, m, n, fptr, frow, fval, fcontrol, &
         finfo, fperm, scale=fscale)
   else
      call f_mc64_matching(job, matrix_type, m, n, fptr, frow, fval, fcontrol, &
         finfo, fperm)
   endif

   ! Copy information out to C structure
   call copy_info_out(finfo, cinfo)

   ! Fix output arrays
   if(.not.f_arrays) then
      call C_F_POINTER(cperm, fperm, shape = (/ m+n /))
      do i = 1, m+n
         if(fperm_alloc(i) < 0) then
            fperm(i) = fperm_alloc(i) + 1
         else
            fperm(i) = fperm_alloc(i) - 1
         endif
      end do
   endif
end subroutine mc64_matching_d
